home *** CD-ROM | disk | FTP | other *** search
/ Graphics Plus / Graphics Plus.iso / general / procssng / ccs / ccs-11tl.lha / lbl / hips / sources / isobuild / iso_march.c < prev    next >
Encoding:
C/C++ Source or Header  |  1992-02-11  |  18.6 KB  |  735 lines
  1. /* iso_march.c    routine for computing the marching cubes triangles
  2.           and edge list */
  3.  
  4. /*  for use with the isobuild program
  5.  *
  6.  * by Brian Tierney, LBL  12/90
  7.  *            Lawrence Berkeley Laboratory
  8.  *            Imaging Technologies Group
  9.  *            email: bltierney@lbl.gov
  10.  *
  11. */
  12.  
  13. /* $Id: iso_march.c,v 1.6 1992/01/31 02:05:45 tierney Exp $ */
  14.  
  15. /* $Log: iso_march.c,v $
  16.  * Revision 1.6  1992/01/31  02:05:45  tierney
  17.  * y
  18.  *
  19.  * Revision 1.5  1992/01/30  20:05:03  davidr
  20.  * prior to Brian's changes
  21.  *
  22.  * Revision 1.4  1992/01/10  01:59:31  davidr
  23.  * works with triserv now
  24.  *
  25.  * Revision 1.2  1991/12/19  01:42:20  davidr
  26.  * added RCS identification markers
  27.  * */
  28.  
  29. static char rcsid[] = "$Id: iso_march.c,v 1.6 1992/01/31 02:05:45 tierney Exp $";
  30.  
  31. #include "isobuild.h"
  32. #include "cell_table.h"
  33.  
  34. /************ Temporary Globals *****************/
  35.  
  36. float     DATA1, DATA2, DATA3, DATA4, DATA5, DATA6, DATA7, DATA8;
  37.  
  38. /**************************** iso_surface ****************************/
  39.  
  40. int
  41. iso_march(startx, starty, endx, endy, slice)
  42.     int       startx, starty, endx, endy, slice;
  43. {
  44.     register int x, y, xdim1, ydim1;
  45.     int       index;
  46.     int       npolys = 0;
  47.     CUBE_TRIANGLES **tmpslice;
  48.     float     crossings[13][3];
  49.  
  50.     for (y = 0; y < 3; y++)    /* initialize triangles array */
  51.     for (x = 0; x < 13; x++)
  52.         crossings[x][y] = 0.;
  53.  
  54.     xdim1 = endx - 1;
  55.     ydim1 = endy - 1;
  56.  
  57.     for (y = starty; y < ydim1; y++) {
  58.     for (x = startx; x < xdim1; x++) {
  59.         index = calc_index_and_temps(x, y, slice);
  60.  
  61.         if (index > 0) {    /* index of 0 means no triangles */
  62.  
  63.         if ((get_cell_verts(index, (float) x, (float) y,
  64.                     (float) slice, (float) TVAL,
  65.                     crossings)) < 0)
  66.             return (-1);
  67.         npolys += get_cell_polys(index, crossings, x, y, slice);
  68.         }
  69.     }
  70.     }
  71.     if (DUP_CHECK) {
  72.     tmpslice = currslice;
  73.     currslice = prevslice;
  74.     prevslice = tmpslice;
  75.     }
  76.     if (VERBOSE2)
  77.     fprintf(stderr, "Slice %d: %d polygons \n", slice, npolys);
  78.  
  79.     return (npolys);
  80. }
  81.  
  82. /********************************************************/
  83. int
  84. iso_march_block(xloc, yloc, zloc, width, height,
  85.         sx, sy, ex, ey)
  86.     int       xloc, yloc, zloc, width, height;
  87.     int       sx, sy, ex, ey;
  88. {
  89.     register int i, j, x, y;
  90.     int       index, npolys = 0;
  91.     float     crossings[13][3];
  92.     CUBE_TRIANGLES **tmpslice;
  93.  
  94.     for (y = 0; y < 3; y++)    /* initialize triangles array */
  95.     for (x = 0; x < 13; x++)
  96.         crossings[x][y] = 0.;
  97.  
  98.     y = yloc;
  99.     for (i = 0; i < height; i++) {
  100.     x = xloc;
  101.     for (j = 0; j < width; j++) {
  102.  
  103.         if (x >= sx && y >= sy && x < ex && y < ey) {
  104.  
  105.         index = calc_index_and_temps(x, y, zloc);
  106.  
  107.         if (index > 0) {/* index of 0 means no triangles */
  108.             if (DUP_CHECK) {    /* initialize */
  109.             currslice[y][x].edge_index = -1;
  110.             currslice[y][x].num_edges = 0;
  111.             }
  112.             if ((get_cell_verts(index, (float) x, (float) y,
  113.                     (float) zloc, (float) TVAL,
  114.                     crossings)) < 0)
  115.             return (-1);
  116.             npolys += get_cell_polys(index, crossings, x, y, zloc);
  117.         }
  118.         }
  119.         x++;
  120.     }
  121.     y++;
  122.     }
  123.     if (DUP_CHECK) {
  124.     tmpslice = currslice;
  125.     currslice = prevslice;
  126.     prevslice = tmpslice;
  127.     }
  128.     return (npolys);
  129. }
  130.  
  131. /************************* calc_index_and_temps ****************************/
  132. int
  133. calc_index_and_temps(x1, y1, z1)
  134.     int       x1, y1, z1;
  135.  
  136. /* This subroutine calculates the index into the cell_table,
  137.  * and creates some global temporary variables (for speed).
  138.  */
  139. {
  140.     int       x2, y2, z2, index = 0;
  141.     Grid_type grid_cube[8];
  142.     float     nval;
  143.  
  144.     x2 = x1 + 1;
  145.     y2 = y1 + 1;
  146.     z2 = z1 + 1;
  147.  
  148.     grid_cube[0] = grid[z1][y1][x1];
  149.     grid_cube[1] = grid[z1][y1][x2];
  150.     grid_cube[2] = grid[z1][y2][x2];
  151.     grid_cube[3] = grid[z1][y2][x1];
  152.     grid_cube[4] = grid[z2][y1][x1];
  153.     grid_cube[5] = grid[z2][y1][x2];
  154.     grid_cube[6] = grid[z2][y2][x2];
  155.     grid_cube[7] = grid[z2][y2][x1];
  156.  
  157.  
  158.     if (grid_cube[0] == 1)
  159.     index++;
  160.     if (grid_cube[1] == 1)
  161.     index += 2;
  162.     if (grid_cube[2] == 1)
  163.     index += 4;
  164.     if (grid_cube[3] == 1)
  165.     index += 8;
  166.  
  167.     if (grid_cube[4] == 1)
  168.     index += 16;
  169.     if (grid_cube[5] == 1)
  170.     index += 32;
  171.     if (grid_cube[6] == 1)
  172.     index += 64;
  173.     if (grid_cube[7] == 1)
  174.     index += 128;
  175.  
  176.     /* if index = 0 or 255, then the surface does not intersect this voxel */
  177.     if (index == 0 || index == 255)
  178.     return (0);
  179.  
  180.     DATA1 = (float) data[z1][y1][x1];
  181.     DATA2 = (float) data[z1][y1][x2];
  182.     DATA3 = (float) data[z1][y2][x2];
  183.     DATA4 = (float) data[z1][y2][x1];
  184.  
  185.     DATA5 = (float) data[z2][y1][x1];
  186.     DATA6 = (float) data[z2][y1][x2];
  187.     DATA7 = (float) data[z2][y2][x2];
  188.     DATA8 = (float) data[z2][y2][x1];
  189.  
  190.  
  191.  
  192.     if (SEG_METHOD > 0 || SMALLER_BOX == 1 || SMOOTH_GRID == 1) {
  193.     /* create false value N % from the thresh value */
  194.     nval = (float) (TVAL - ((data_max - data_min) * .1));
  195.  
  196.     if (DATA1 > TVAL && grid_cube[0] <= 0)
  197.         DATA1 = nval;
  198.  
  199.     if (DATA2 > TVAL && grid_cube[1] <= 0)
  200.         DATA2 = nval;
  201.  
  202.     if (DATA3 > TVAL && grid_cube[2] <= 0)
  203.         DATA3 = nval;
  204.  
  205.     if (DATA4 > TVAL && grid_cube[3] <= 0)
  206.         DATA4 = nval;
  207.  
  208.     if (DATA5 > TVAL && grid_cube[4] <= 0)
  209.         DATA5 = nval;
  210.  
  211.     if (DATA6 > TVAL && grid_cube[5] <= 0)
  212.         DATA6 = nval;
  213.  
  214.     if (DATA7 > TVAL && grid_cube[6] <= 0)
  215.         DATA7 = nval;
  216.  
  217.     if (DATA8 > TVAL && grid_cube[7] <= 0)
  218.         DATA8 = nval;
  219.     }
  220.  
  221.     return (index);
  222. }
  223.  
  224. /**************************** get_cell_verts ****************************/
  225.  
  226. int
  227. get_cell_verts(index, x1, y1, z1, threshold, crossings)
  228.     int       index;
  229.     float     x1, y1, z1;
  230.     float     threshold;
  231.     float     crossings[13][3];
  232. {
  233. /* This routine computes the vertex locations */
  234.  
  235.     register int i;
  236.     float     x2, y2, z2, val;
  237.     int       nedges;
  238.     int       crnt_edge;
  239.  
  240. #define linterp(d1,d2,t,x) ((t-d1) / ((d2-d1)) + x)
  241.  
  242.     x2 = x1 + 1.;
  243.     y2 = y1 + 1.;
  244.     z2 = z1 + 1.;
  245.  
  246.     nedges = cell_table[index].nedges;
  247.     for (i = 0; i < nedges; i++) {
  248.     crnt_edge = cell_table[index].edges[i];
  249.     switch (crnt_edge) {
  250.     case 1:
  251.         if (DATA2 == DATA1)
  252.         val = x1;
  253.         else {
  254.         val = linterp(DATA1, DATA2, threshold, x1);
  255.         }
  256.         crossings[1][0] = val;
  257.         crossings[1][1] = y1;
  258.         crossings[1][2] = z1;
  259.         break;
  260.     case 2:
  261.         if (DATA3 == DATA2)
  262.         val = y1;
  263.         else {
  264.         val = linterp(DATA2, DATA3, threshold, y1);
  265.         }
  266.         crossings[2][1] = val;
  267.         crossings[2][0] = x2;
  268.         crossings[2][2] = z1;
  269.         break;
  270.  
  271.     case 3:
  272.         if (DATA3 == DATA4)
  273.         val = x1;
  274.         else {
  275.         val = linterp(DATA4, DATA3, threshold, x1);
  276.         }
  277.         crossings[3][0] = val;
  278.         crossings[3][1] = y2;
  279.         crossings[3][2] = z1;
  280.         break;
  281.  
  282.     case 4:
  283.         if (DATA4 == DATA1)
  284.         val = y1;
  285.         else {
  286.         val = linterp(DATA1, DATA4, threshold, y1);
  287.         }
  288.         crossings[4][1] = val;
  289.         crossings[4][0] = x1;
  290.         crossings[4][2] = z1;
  291.         break;
  292.     case 5:
  293.         if (DATA6 == DATA5)
  294.         val = x1;
  295.         else {
  296.         val = linterp(DATA5, DATA6, threshold, x1);
  297.         }
  298.         crossings[5][0] = val;
  299.         crossings[5][1] = y1;
  300.         crossings[5][2] = z2;
  301.         break;
  302.  
  303.     case 6:
  304.         if (DATA7 == DATA6)
  305.         val = y1;
  306.         else {
  307.         val = linterp(DATA6, DATA7, threshold, y1);
  308.         }
  309.         crossings[6][1] = val;
  310.         crossings[6][0] = x2;
  311.         crossings[6][2] = z2;
  312.         break;
  313.  
  314.     case 7:
  315.         if (DATA7 == DATA8)
  316.         val = x1;
  317.         else {
  318.         val = linterp(DATA8, DATA7, threshold, x1);
  319.         }
  320.         crossings[7][0] = val;
  321.         crossings[7][1] = y2;
  322.         crossings[7][2] = z2;
  323.         break;
  324.  
  325.     case 8:
  326.         if (DATA8 == DATA5)
  327.         val = y1;
  328.         else {
  329.         val = linterp(DATA5, DATA8, threshold, y1);
  330.         }
  331.         crossings[8][1] = val;
  332.         crossings[8][0] = x1;
  333.         crossings[8][2] = z2;
  334.         break;
  335.  
  336.     case 9:
  337.         if (DATA5 == DATA1)
  338.         val = z1;
  339.         else {
  340.         val = linterp(DATA1, DATA5, threshold, z1);
  341.         }
  342.         crossings[9][2] = val;
  343.         crossings[9][1] = y1;
  344.         crossings[9][0] = x1;
  345.         break;
  346.  
  347.     case 10:
  348.         if (DATA6 == DATA2)
  349.         val = z1;
  350.         else {
  351.         val = linterp(DATA2, DATA6, threshold, z1);
  352.         }
  353.         crossings[10][2] = val;
  354.         crossings[10][1] = y1;
  355.         crossings[10][0] = x2;
  356.         break;
  357.  
  358.     case 11:
  359.         if (DATA8 == DATA4)
  360.         val = z1;
  361.         else {
  362.         val = linterp(DATA4, DATA8, threshold, z1);
  363.         }
  364.         crossings[11][2] = val;
  365.         crossings[11][1] = y2;
  366.         crossings[11][0] = x1;
  367.         break;
  368.  
  369.     case 12:
  370.         if (DATA7 == DATA3)
  371.         val = z1;
  372.         else {
  373.         val = linterp(DATA3, DATA7, threshold, z1);
  374.         }
  375.         crossings[12][2] = val;
  376.         crossings[12][1] = y2;
  377.         crossings[12][0] = x2;
  378.         break;
  379.  
  380.     }            /* end switch */
  381.     }                /* end for */
  382.  
  383. #define TEST_FOR_NEG_VALUES
  384.  
  385. #ifdef TEST_FOR_NEG_VALUES
  386.     for (i = 0; i < 13; i++) {
  387.     if (crossings[i][0] < 0) {
  388.         fprintf(stderr,
  389.             "Neg val(%f) at crossing[%d][0], from loc %d,%d,%d, index=%d, thresh=%d \n",
  390.             crossings[i][0], i, (int) x1, (int) y1, (int) z1, index, (int) threshold);
  391.         fprintf(stderr, " DATA values: %f,%f,%f,%f,%f,%f,%f,%f \n",
  392.             DATA1, DATA2, DATA3, DATA4, DATA5, DATA6, DATA7, DATA8);
  393.         /* crossings[i][0] = 0.; */
  394.     }
  395.     if (crossings[i][1] < 0) {
  396.         fprintf(stderr,
  397.             "Neg val(%f) at crossing[%d][0], from loc %d,%d,%d, index=%d, thresh=%d \n",
  398.             crossings[i][1], i, (int) x1, (int) y1, (int) z1, index, (int) threshold);
  399.         fprintf(stderr, " DATA values: %f,%f,%f,%f,%f,%f,%f,%f \n",
  400.             DATA1, DATA2, DATA3, DATA4, DATA5, DATA6, DATA7, DATA8);
  401.         /* crossings[i][1] = 0.; */
  402.     }
  403.     if (crossings[i][2] < 0) {
  404.         fprintf(stderr,
  405.             "Neg val(%f) at crossing[%d][0], from loc %d,%d,%d, index=%d, thresh=%d \n",
  406.             crossings[i][2], i, (int) x1, (int) y1, (int) z1, index, (int) threshold);
  407.         fprintf(stderr, " DATA values: %f,%f,%f,%f,%f,%f,%f,%f \n",
  408.             DATA1, DATA2, DATA3, DATA4, DATA5, DATA6, DATA7, DATA8);
  409.         /* crossings[i][2] = 0.; */
  410.     }
  411.     }
  412. #endif
  413.  
  414.     for (i = 0; i < 13; i++)
  415.     if (crossings[i][0] < 0 || crossings[i][1] < 0 || crossings[i][2] < 0) {
  416.         Error("negative vertex detected");
  417.     }
  418.     return (0);
  419. }
  420.  
  421. /**************************** get_cell_polys ****************************/
  422.  
  423. int
  424. get_cell_polys(index, crossings, x, y, z)
  425.     int       index;
  426.     float     crossings[13][3];
  427.     int       x, y, z;
  428. /* This subroutine will calculate the polygons */
  429. {
  430.     register int num_polys, poly_cnt;
  431.     register int poly, idx;
  432.     float    *p1, *p2, *p3;
  433.     void      add_polygon();
  434.  
  435.     num_polys = cell_table[index].npolys;    /* possible values: 0 to 10 */
  436.  
  437.     idx = poly_cnt = 0;
  438.     for (poly = 0; poly < num_polys; poly++) {
  439.  
  440.     p1 = &crossings[cell_table[index].polys[idx++]][0];
  441.     p2 = &crossings[cell_table[index].polys[idx++]][0];
  442.     p3 = &crossings[cell_table[index].polys[idx++]][0];
  443.  
  444. /* #define SHOW_TRIANGLES  */
  445. #ifdef SHOW_TRIANGLES
  446.     fprintf(stderr, "triangle from data cube at (%d,%d,%d): \n", x, y, z);
  447.     fprintf(stderr, "  (%.3f,%.3f,%.3f) (%.3f,%.3f,%.3f) (%.3f, %.3f,%.3f) \n",
  448.          p1[0], p1[1], p1[2], p2[0], p2[1], p2[2], p3[0], p3[1], p3[2]);
  449. #endif
  450.  
  451. #define CHECK_FOR_DEGENERATE_TRIANGLES
  452.  
  453. #ifdef CHECK_FOR_DEGENERATE_TRIANGLES
  454.     /*
  455.      * degenerate triangles are produced when one of the corners of the
  456.      * cube equals the threshold value, and some renderers don't handle
  457.      * degenerate triangles well
  458.      */
  459.  
  460.     if ((p1[0] == p2[0] && p1[1] == p2[1] && p1[2] == p2[2]) ||
  461.         (p1[0] == p3[0] && p1[1] == p3[1] && p1[2] == p3[2]) ||
  462.         (p2[0] == p3[0] && p2[1] == p3[1] && p2[2] == p3[2]));
  463.     else {
  464. #endif
  465.         /* p1,p2,p3 contain x,y,z values for each vertex of the triangle */
  466.         add_polygon(p1, p2, p3, x, y, z);
  467.         poly_cnt++;
  468. #ifdef CHECK_FOR_DEGENERATE_TRIANGLES
  469.     }
  470. #endif
  471.     }
  472.  
  473.     return (poly_cnt);
  474. }
  475.  
  476. /********************************************************************/
  477.  
  478. /* if SMALL_CHUNKS is defined, output will be written whenever the
  479.    buffer is full, otherwise the size of the buffer is increased
  480.  */
  481.  
  482. VERT_PTR  vertex_list = NULL;
  483. NORM_PTR  norm_list = NULL;
  484. int      *conn_list = NULL;
  485.  
  486. #define OUTBUF_SIZE 50000
  487.  
  488. int       VERT_LIMIT = OUTBUF_SIZE;
  489. int       CONN_LIMIT = OUTBUF_SIZE * 3;
  490.  
  491. #define VERT_INCR 30000        /* size to increase vertex lists when full */
  492.  
  493. static int num_vertices = 0, num_conn = 0;
  494.  
  495. /**************************** add_polygon ****************************/
  496.  
  497. void
  498. add_polygon(p1, p2, p3, x, y, z)
  499.     float    *p1, *p2, *p3;
  500.     int       x, y, z;
  501. {
  502.     int       idx, conn1, conn2, conn3;
  503.     NORMAL_VECT norm, calc_normal();
  504.     int       size;
  505.  
  506.     if (vertex_list == NULL)
  507.     vertex_list = (VERT_PTR) malloc(sizeof(VERTEX) * VERT_LIMIT);
  508.  
  509.     if (MARCH_NORMALS) {
  510.     if (norm_list == NULL)
  511.         norm_list = (NORM_PTR) malloc(sizeof(NORMAL) * VERT_LIMIT);
  512.     }
  513.     if (DUP_CHECK && conn_list == NULL) {
  514.     conn_list = (int *) malloc(sizeof(int) * CONN_LIMIT);
  515.     }
  516.     /* store the vertices */
  517.     if (!DUP_CHECK || (idx = check_for_duplicate_vertex(p1, x, y, z)) < 0) {
  518.     vertex_list[num_vertices].x = p1[0];
  519.     vertex_list[num_vertices].y = p1[1];
  520.     vertex_list[num_vertices].z = p1[2];
  521.  
  522.     if (MARCH_NORMALS) {
  523.         if (PRE_NORMALS) {
  524.         norm_list[num_vertices].nx = normals[z][y][x].x;
  525.         norm_list[num_vertices].ny = normals[z][y][x].y;
  526.         norm_list[num_vertices].nz = normals[z][y][x].z;
  527.         } else {
  528.         norm = calc_normal((int) (p1[0] + .5), (int) (p1[1] + .5),
  529.                    (int) (p1[2] + .5));
  530.         norm_list[num_vertices].nx = norm.x;
  531.         norm_list[num_vertices].ny = norm.y;
  532.         norm_list[num_vertices].nz = norm.z;
  533.         }
  534.     }
  535.     conn1 = num_vertices++;
  536.     } else {
  537.     conn1 = idx;
  538.     }
  539.  
  540.     if (!DUP_CHECK || (idx = check_for_duplicate_vertex(p2, x, y, z)) < 0) {
  541.     vertex_list[num_vertices].x = p2[0];
  542.     vertex_list[num_vertices].y = p2[1];
  543.     vertex_list[num_vertices].z = p2[2];
  544.  
  545.     if (MARCH_NORMALS) {
  546.         if (PRE_NORMALS) {
  547.         norm_list[num_vertices].nx = normals[z][y][x].x;
  548.         norm_list[num_vertices].ny = normals[z][y][x].y;
  549.         norm_list[num_vertices].nz = normals[z][y][x].z;
  550.         } else {
  551.         norm = calc_normal((int) (p2[0] + .5), (int) (p2[1] + .5),
  552.                    (int) (p2[2] + .5));
  553.         norm_list[num_vertices].nx = norm.x;
  554.         norm_list[num_vertices].ny = norm.y;
  555.         norm_list[num_vertices].nz = norm.z;
  556.         }
  557.     }
  558.     conn2 = num_vertices++;
  559.     } else {
  560.     conn2 = idx;
  561.     }
  562.  
  563.     if (!DUP_CHECK || (idx = check_for_duplicate_vertex(p3, x, y, z)) < 0) {
  564.     vertex_list[num_vertices].x = p3[0];
  565.     vertex_list[num_vertices].y = p3[1];
  566.     vertex_list[num_vertices].z = p3[2];
  567.  
  568.     if (MARCH_NORMALS) {
  569.         if (PRE_NORMALS) {
  570.         norm_list[num_vertices].nx = normals[z][y][x].x;
  571.         norm_list[num_vertices].ny = normals[z][y][x].y;
  572.         norm_list[num_vertices].nz = normals[z][y][x].z;
  573.         } else {
  574.         norm = calc_normal((int) (p3[0] + .5), (int) (p3[1] + .5),
  575.                    (int) (p3[2] + .5));
  576.         norm_list[num_vertices].nx = norm.x;
  577.         norm_list[num_vertices].ny = norm.y;
  578.         norm_list[num_vertices].nz = norm.z;
  579.         }
  580.     }
  581.     conn3 = num_vertices++;
  582.     } else {
  583.     conn3 = idx;
  584.     }
  585.  
  586.     if (DUP_CHECK) {
  587.     /* store the connnectity info info */
  588.     if (currslice[y][x].edge_index == -1)
  589.         currslice[y][x].edge_index = num_conn;
  590.     currslice[y][x].num_edges += 3;
  591.  
  592.     conn_list[num_conn++] = conn1;
  593.     conn_list[num_conn++] = conn2;
  594.     conn_list[num_conn++] = conn3;
  595.     }
  596.     if (SMALL_CHUNKS) {
  597.     if (num_vertices + 4 >= VERT_LIMIT || num_conn + 3 >= CONN_LIMIT) {
  598.         dump_polys(1);    /* write polys found so far */
  599.  
  600.         if (DUP_CHECK)
  601.         init_dup_vertex_check_arrays();
  602.     }
  603.     } else {
  604.     if (num_vertices + 4 >= VERT_LIMIT) {
  605.         /* get more space */
  606.         VERT_LIMIT += VERT_INCR;    /* add this much space */
  607.         size = VERT_LIMIT * sizeof(VERTEX);    /* size for malloc/realloc */
  608.         if ((vertex_list = (VERT_PTR) realloc((char *) vertex_list,
  609.                           size)) == NULL) {
  610.         Error("not enough memory to store vertices");
  611.         }
  612.         if (MARCH_NORMALS) {
  613.         size = VERT_LIMIT * sizeof(NORMAL);    /* size for
  614.                              * malloc/realloc */
  615.         if ((norm_list = (NORM_PTR) realloc((char *) norm_list,
  616.                             size)) == NULL) {
  617.             Error("not enough memory to store normals");
  618.         }
  619.         }
  620.     }
  621.     if (DUP_CHECK && num_conn + 3 >= CONN_LIMIT) {
  622.         CONN_LIMIT += VERT_INCR * 3;    /* add this much space */
  623.         size = CONN_LIMIT * sizeof(int);
  624.         if ((conn_list = (int *) realloc((char *) conn_list,
  625.                          size)) == NULL) {
  626.         Error("not enough memory to store connectivity");
  627.         }
  628.     }
  629.     }
  630. }
  631.  
  632. /***********  check for duplicate vertices  *******************************/
  633. int
  634. check_for_duplicate_vertex(vert, x, y, z)
  635.     float     vert[3];        /* vertex location */
  636.     int       x, y, z;        /* location in original data set */
  637. {
  638.     int       idx;
  639.  
  640.  
  641.     /* first, check for duplicate vertices from the same cube */
  642.     if ((idx = search_for_dup(vert, currslice[y][x].num_edges,
  643.                   currslice[y][x].edge_index)) >= 0)
  644.     return (idx);
  645.  
  646.     if (x > 0) {
  647.     if ((idx = search_for_dup(vert, currslice[y][x - 1].num_edges,
  648.                   currslice[y][x - 1].edge_index)) >= 0)
  649.         return (idx);
  650.     }
  651.     if (y > 0) {
  652.     if ((idx = search_for_dup(vert, currslice[y - 1][x].num_edges,
  653.                   currslice[y - 1][x].edge_index)) >= 0)
  654.         return (idx);
  655.     }
  656.     if (z > 0) {
  657.     if ((idx = search_for_dup(vert, prevslice[y][x].num_edges,
  658.                   prevslice[y][x].edge_index)) >= 0)
  659.         return (idx);
  660.     }
  661.     return (-1);
  662. }
  663.  
  664. /***************************************************************/
  665. int
  666. search_for_dup(vert, num_to_check, conn_idx)
  667.     float     vert[3];        /* vertex location */
  668.     int       num_to_check;    /* number of vertices in this cube */
  669.     int       conn_idx;        /* index into connectivity list */
  670. {
  671.     register int i, idx;
  672.     static int dup_cnt = 0;
  673.  
  674.     if (conn_idx == -1)
  675.     return (-1);
  676.  
  677.     if (num_to_check == -1) {    /* for testing */
  678.     fprintf(stderr, "Found %d duplicate vertices \n", dup_cnt);
  679.     dup_cnt = 0;
  680.     return (-1);
  681.     }
  682.     for (i = 0; i < num_to_check; i++) {
  683.     idx = conn_list[conn_idx];
  684.     if (idx > num_vertices) {    /* seems to need this check, but I
  685.                      * dont know why? */
  686.         fprintf(stderr, "Error checking dup vertex: conn_idx = %d, num_verts = %d \n",
  687.             conn_idx, num_vertices);
  688.         return (-1);
  689.     }
  690.     if (vertex_list[idx].x == vert[0] && vertex_list[idx].y == vert[1] &&
  691.         vertex_list[idx].z == vert[2]) {
  692.         dup_cnt++;
  693.         return (idx);
  694.     }
  695.     conn_idx++;
  696.     }
  697.     return (-1);        /* none found */
  698. }
  699.  
  700. /**************************** dump_polys ****************************/
  701. int
  702. dump_polys(more)
  703.     int       more;
  704. {
  705.     void      write_polys_to_socket();
  706.  
  707.     if (VERBOSE2)
  708.     Status("Writing triangles to file/socket..");
  709.  
  710.     if (SERVER) {
  711.     write_polys_to_socket(vertex_list, norm_list, conn_list,
  712.                   num_vertices, num_conn, more);
  713.     } else if (OUTPUT_TYPE == 0) {
  714.     write_polys(vertex_list, norm_list, conn_list, num_vertices,
  715.             num_conn, more);
  716.     } else if (OUTPUT_TYPE == 1) {
  717.     write_polys_ascii(vertex_list, norm_list, conn_list,
  718.               num_vertices, num_conn);
  719.     } else if (OUTPUT_TYPE == 2) {
  720.     write_polys_movie_byu(vertex_list, conn_list, num_vertices, num_conn);
  721.     } else if (OUTPUT_TYPE == 3) {
  722.     write_polys_wavefront(vertex_list, norm_list, conn_list,
  723.                   num_vertices, num_conn);
  724.     } else if (OUTPUT_TYPE == 4) {
  725.     write_polys_sunvision(vertex_list, norm_list, conn_list,
  726.                   num_vertices, num_conn);
  727.     }
  728.     num_vertices = num_conn = 0;
  729.  
  730.     if (VERBOSE)
  731.     search_for_dup(0, -1, 0);    /* prints diag message */
  732.  
  733.     return;
  734. }
  735.